library(tidyverse)
# library(rvest)
# library(albersusa)
library(sf)
library(leaflet)
library(leafpop)
library(formattable)
library(ggiraph)
library(urbnmapr) # change to include AK/HI on map

# library(covdata)

# Default ggplot theme
theme_set(theme_bw() +
            theme(axis.title.x = element_text(size = rel(1.25)),
                  axis.title.y = element_text(size = rel(1.25)),
                  plot.caption = element_text(face = "italic", size = 9),
                  legend.position = "bottom",
                  legend.background = element_rect(fill = "transparent"), 
                  panel.background = element_rect(fill = "transparent"),
                  panel.border = element_blank(),
                  axis.line = element_line(color="black"),
                  panel.grid.major.x = element_blank(),
                  panel.grid.minor.x = element_blank(),
                  panel.grid.minor.y = element_blank(),
                  panel.grid.major.y = element_line(color = "grey90", linetype = 2)
            ))

Last updated: 2021-07-06 14:07:39 All source documents for this report can be found on github

VACCINE SIGN-UP

If you are 12 years old or older, you are eligible to get your vaccine in Nashville (for free)- Go to one of these sites to sign up:

Disclaimer

Update 5/10/2021 New Cases now reported only by the date of detection. Will update with details once it becomes clear how these results will be mixed with the prior reporting method (date of report).

The Tennessee state information is simply a visual summary of the data directly downloaded from TN state Dept of Health, by date reported (rather than date of test). I do not attempt to model, make predictions, or give advice on how to use this data. If you want that sort of thing, go to the VUMC Health Policy site, which has developed some very good local models and Advisory Memos for Tennessee.

For the most up-to-date numbers visit the TN state Dept of Health or local health department sites. In particular, the Metro. Nashville site has a very nice summary- and note that their report of daily cases usually precedes the State Health Dept, due to the verification process and assignment of county of residence at the state.

Other useful sites for following trends in Tennessee and nationally:

TN Dept of Health Data

These data are taken directly from the TN state Dept of Health. Previously the DOH reported results by date of report, rather than date of test. They transitioned to reporting results by Date of Test on 5/10/2021, although it is not clear if data reported prior to that will be cleaned up.

# READ IN AND format data ####
tndoh_daily <- readxl::read_excel(here::here("data_tndoh/Public-Dataset-Daily-Case-Info.xlsx"))
names(tndoh_daily) <- tolower(names(tndoh_daily))

# remove the new_test results which our bulk reported or adjusted
tndoh_daily <-
  tndoh_daily %>%  
  mutate(date = lubridate::ymd(date)) %>%  
  arrange(date) 

# tndoh_daily <-
#   tndoh_daily %>%
#   filter(date > lubridate::ymd("2020-01-01")) # data error  introduced; filter out erroneous entry for 1960-01-01
   
tndoh_daily[tndoh_daily$date %in% lubridate::ymd(c("2020-07-13", "2020-07-11", "2020-06-12")), "new_tests"] <- NA_real_

# STATE POPULATION DATA ####
# us_data <-  readRDS(here::here("data_rds/us_data.rds"))
us_data <- get_urbn_map("states", sf = TRUE)
us_data <- 
  us_data %>% 
  mutate(name = state_name,
         stusps = state_abbv)


# state_pop <- readxl::read_excel(here::here("state_populations.xlsx"))
# state_pop$region <- gsub("\\.", "", state_pop$region)
# state_pop <-
#   state_pop %>% 
#   pivot_longer(cols = -c(region, census_2010),
#                names_to = "pop_year",
#                names_prefix = "pop_",
#                values_to = "pop_estimate")
# saveRDS(state_pop, here::here("data_rds/state_pop.rds"))

state_pop <- readRDS( here::here("data_rds/state_pop.rds"))
state_pop_2019 <-
  state_pop %>%  
  filter(pop_year == 2019) %>%  
  filter(region %in% state.name)
TN_pop2019 <- state_pop_2019$pop_estimate[state_pop_2019$region == "Tennessee"]

# vaccine state info
tndoh_vax_state <- readxl::read_excel(here::here("data_tndoh/COVID_VACCINE_STATE_SUMMARY.XLSX"))
names(tndoh_vax_state) <- tolower(names(tndoh_vax_state))
tndoh_vax_state <-
  tndoh_vax_state %>% 
  mutate(date = lubridate::ymd(date),
         vaccine_count_pct = vaccine_count/TN_pop2019*100,
         recipient_count_pct = recipient_count/TN_pop2019*100,
         recip_fully_vacc_pct = recip_fully_vacc/TN_pop2019*100) 

# create annotations
tndoh_daily <-
  left_join(tndoh_daily, tndoh_vax_state,
            by = "date") %>%  
  mutate(pos_tests_new = pos_tests -lag(pos_tests, 1), 
         pos_pct =   pos_tests/total_tests,
         pos_pct_today = pos_tests_new/new_tests,
         hospitalized_tot = replace_na(total_hosp, 0),
         hospitalized_new = total_hosp - lag(total_hosp),
         new_cases_100k = new_cases/TN_pop2019*100000,
         new_cases_7davg = RcppRoll::roll_meanr(new_cases, 7),
         new_cases_7davg_100k = new_cases_7davg/TN_pop2019*100000,
         new_deaths_7davg = RcppRoll::roll_meanr(new_deaths, 7),
         new_hosp_7davg = RcppRoll::roll_meanr(new_hosp, 7),
         new_tests_7davg = RcppRoll::roll_meanr(new_tests, 7),
         pos_pct_7davg = RcppRoll::roll_meanr(pos_tests_new, 7)/RcppRoll::roll_meanr(new_tests, 7),
         tt_daily = sprintf("<tr><td><b>%s</b></td></tr><br>
                     <tr><td><b>New Cases: %s</b></td></tr><br>
                     <tr><td><b>New Cases/100k: %s</b></td></tr><br>
                     <tr><td><b>New Hospitalizations: %s</b></td></tr><br>
                     <tr><td><b>New Deaths: %s</b></td></tr><br>
                     <tr><td><b>Test+ %% Today: %1.1f%%</b></td></tr><br>
                     <tr><td>New Cases 7d Avg/100k: %s</td></tr><br>
                     <tr><td>Deaths 7d Avg: %1.1f</td></tr><br>
                     <tr><td>Test+ %% 7d Avg: %1.1f%%</td></tr><br>
                     <tr><td>Case Fatality Ratio: %1.2f%%</td></tr><br>
                     <tr><td>New Vaccinations: %s</td></tr><br>
                     <tr><td>Rec 1st Vaccination %%Pop: %1.2f%%</td></tr>",
                     paste0(month.name[lubridate::month(date)], " ", lubridate::mday(date), ", ", lubridate::year(date), " (", lubridate::wday(date, label = T), ")"), 
                     comma(new_cases, digits = 0), 
                     comma(new_cases_100k, digits = 0), 
                     comma(new_hosp, 0), 
                     comma(new_deaths, digits = 0), 
                     pos_pct_today*100, 
                     comma(new_cases_7davg_100k, 1),
                     comma(new_deaths_7davg, digits = 0), 
                     pos_pct_7davg*100, 
                     total_deaths/total_cases*100,
                     comma(new_vaccine_count, 0),
                     recipient_count_pct) ) 


tndoh_daily_limited <-
  tndoh_daily %>% 
  select(date, total_cases, total_hosp, total_deaths) %>%  
  pivot_longer(cols = c(total_cases, total_hosp, total_deaths),
               names_to = "measure",
               values_to = "count") %>%  
  mutate(var.f = factor(measure, levels = c("total_cases", "total_hosp", "total_deaths"),
                        labels = c("Cases", "Hospitalizations", "Deaths")))  
title_tnsdc <- cowplot::ggdraw() + 
  cowplot::draw_label(paste0("COVID-19 in Tennessee (as of ", max(tndoh_daily$date), ")") , fontface='bold')

tndoh_daily_long <-
  tndoh_daily %>%  
  select(date, new_cases, new_deaths, new_hosp) %>%   
  pivot_longer(cols = c(new_cases, new_deaths, new_hosp),
               names_to = "var", 
               values_to = "value") %>%  
  mutate(var.f = factor(var, levels = c("new_cases", "new_deaths", "new_hosp"), 
                        labels = c("Cases", "Hospitalizations", "Deaths")))

label_pos_pct <-
    with(tndoh_daily,
         sprintf("Overall Positive: %1.1f%%\nToday's Positive: %1.1f%%", 
                 pos_pct[which(date == max(date))]*100,
                 pos_pct_today[which(date == max(date))]*100))



# AGE DATA ####
tndoh_age <- readxl::read_excel(here::here("data_tndoh/Public-Dataset-Age.xlsx"))
names(tndoh_age) <-  tolower(names(tndoh_age))
tndoh_age <-
  tndoh_age %>% 
  filter(age_range != "Pending") %>%  
  mutate(date = lubridate::ymd(date),
         age_range = gsub(" years", "", age_range),
         age.f = factor(age_range),
         deaths_case_rate = case_when(age.f == "Pending" ~ NA_real_,
                                      TRUE ~ ar_totaldeaths/ar_totalcases*100))

# TN age groups 2019 estimates(taken from https://www.census.gov/data/tables/time-series/demo/popest/2010s-state-detail.html)
tn_age_groups <- 
  tibble::tibble(age.f = factor(c("0-10", "11-20", "21-30", "31-40", "41-50", "51-60", "61-70", "71-80", "81+"),
                                    levels = c("0-10", "11-20", "21-30", "31-40", "41-50", "51-60", "61-70", "71-80", "81+")),
                 populationtn_byage_2019 = c(821984L, 855582L, 942211L, 887898L, 850125L, 893578L, 808164L, 516865L, 252767L)) 

tndoh_age <- 
  left_join(tndoh_age,
          tn_age_groups,
          by = "age.f") %>%  
  group_by(age.f) %>%  
  mutate(ar_new_cases_100k = ar_newcases/populationtn_byage_2019*100000,
         ar_new_cases_100k_7davg = round(RcppRoll::roll_meanr(ar_newcases, 7)/populationtn_byage_2019*100000, 1),
         ar_totaldeaths_100k    = ar_totaldeaths/populationtn_byage_2019*100000,
         ar_new_deaths_100k = ar_newdeaths/populationtn_byage_2019*100000,
         ar_new_deaths_100k_7davg = round(RcppRoll::roll_meanr(ar_newdeaths, 7)/populationtn_byage_2019*100000, 1)) 



# VACCINE by age
 tndoh_vax_age <- readxl::read_excel(here::here("data_tndoh/COVID_VACCINE_AGE_GROUPS.XLSX"),
                                      col_types = c("text", "date", rep("numeric", 6)))
names(tndoh_vax_age) <- tolower(names(tndoh_vax_age))
tndoh_vax_age <-
  tndoh_vax_age %>%
  filter(!grepl("pending", age_group, ignore.case = T )) %>%  
  mutate(date = lubridate::ymd(date),
         age.f = factor(age_group,
                        levels = c("0-10", "11-20", "16-20", "21-30", "31-40", "41-50", "51-60", "61-70", "71-80", "81+"), # 16-20 needs to match other age category
                        labels = c("0-10", "11-20", "11-20", "21-30", "31-40", "41-50", "51-60", "61-70", "71-80", "81+")))
tndoh_age <-
  left_join(tndoh_age, tndoh_vax_age,
            by = c("date", "age.f")) %>%
  rowwise() %>%
  mutate(ar_totalvax_pct = vaccine_count/populationtn_byage_2019*100,
         ar_recipientvax_pct = recipient_count/populationtn_byage_2019*100) %>%
  ungroup() 

# most recent vaccination by age
vax_by_age_mostrecent <-
  tndoh_age %>% 
  filter(!is.na(vaccine_count)) %>% 
  filter(date == max(date, na.rm = T)) %>% 
    select(date, recipient_count, age.f, ar_recipientvax_pct) 
mostrecent_vax_agedate <- max(vax_by_age_mostrecent$date)  

# Deaths by age table ####
age_death_table <-
  tndoh_age %>%
  filter(date == max(date, na.rm = T)) %>% 
  select(age.f,populationtn_byage_2019, ar_totalcases, ar_totaldeaths, ar_totaldeaths_100k) %>% 
  mutate(#ar_totalcases = comma(ar_totalcases, 0),
    ar_totaldeaths_pct = sprintf("%1.2f%%", ar_totaldeaths_100k/1000),
    ar_totaldeaths_100k = sprintf("%1.2f", ar_totaldeaths_100k)) %>% 
  left_join(vax_by_age_mostrecent, by = "age.f")


# COUNTY DATA ####
# tndoh_countynew <- readxl::read_excel(here::here("tn_data_center/TN_COVID19_CountyDaily.xlsx"))
TN_counties <- readRDS(here::here("data_rds/TN_counties.rds")) # see below for source

tndoh_county_coln <-   length(names(readxl::read_excel(here::here("data_tndoh/Public-Dataset-County-New.xlsx"), n_max = 2)))

tndoh_countynew <-  readxl::read_excel(here::here("data_tndoh/Public-Dataset-County-New.xlsx"),
                                       col_types = c("date", "text", rep("numeric", tndoh_county_coln-2))) 

names(tndoh_countynew) <-  tolower(names(tndoh_countynew))
tndoh_countynew <-
  tndoh_countynew %>%  
  mutate(date = lubridate::ymd(date),
         pct_pos = total_cases/(total_tests)*100,
         county_lower = tolower(county)) %>%  
  group_by(county) %>% 
  arrange(date) %>%  
  mutate(new_cases_county_4davg = round(RcppRoll::roll_meanr(new_cases, 4), 1),
         new_cases_county_7davg = round(RcppRoll::roll_meanr(new_cases, 7), 1),
         new_deaths_county_7davg = round(RcppRoll::roll_meanr(new_deaths, 7), 1),
         new_testpos_pct   = round(new_pos_tests/new_tests*100, 1),
         new_testpos_pct_7davg   = round(RcppRoll::roll_meanr(new_pos_tests, 7)/RcppRoll::roll_meanr(new_tests, 7)*100, 1),
         new_hospitalized_county_7davg = round(RcppRoll::roll_meanr(new_hospitalized, 7), 1),
         new_cases_county_7davg_1wkdelta = (new_cases_county_7davg - lag(new_cases_county_7davg, 7))/lag(new_cases_county_7davg, 7) 
  ) %>%  
  ungroup()


  
# TN_data <- tigris::counties(state = "TN", cb = TRUE)
# names(TN_data) <- tolower(names(TN_data))
# TN_counties <- sf::st_as_sf(TN_data) %>%  
#     mutate(county_lower = tolower(as.character(name)))

# TN COUNTY POPULATION FROM US CENSUS 2019 ESTIMATES
# tn_county_pop <- read_csv(here::here("data_csv/TN 2019 county pop estimates.csv"))
# tn_county_pop <-
#   tn_county_pop %>%  
#   mutate(county_lower = tolower(county)) %>%  
#   select(-county)
# TN_counties <-
#   left_join(TN_counties,
#             tn_county_pop,
#             by = "county_lower")


# saveRDS(TN_counties, here::here("data_rds/TN_counties.rds"))
# TN_counties <- readRDS(here::here("data_rds/TN_counties.rds"))

tndoh_countynew <-
  left_join(tndoh_countynew,
            TN_counties,
            by = "county_lower") %>%  
  mutate(total_cases_per100kpop = total_cases/population*100000,
         new_cases_per100kpop_7davg   = new_cases_county_7davg /population*100000,
         total_tests_pctpop = total_tests/population*100,
         new_tests_per100kpop   = new_tests/population*100000,
         total_deaths_per100kpop = total_deaths/population*100000) 

tndoh_countynew_calc <-
  tndoh_countynew %>%  
    group_by(date) %>%  
    filter(tolower(county) %in% TN_counties$county_lower) %>%  # exclude out of state and pending
    mutate(#new_testpos_pct = new_pos_tests/new_tests*100,
           new_cases_county_7davg_rank = rank(-new_cases_county_7davg),
           new_cases_per100kpop_7davg_rank = rank(-new_cases_per100kpop_7davg),
           new_cases_rank = rank(-new_cases)) %>%  
    ungroup() %>%  
    select(county, date, new_cases_county_7davg_rank, 
           new_cases_per100kpop_7davg_rank, new_cases_rank)
tndoh_countynew <-   
  left_join(tndoh_countynew,
            tndoh_countynew_calc,
            by = c("date", "county"))  

 

# DEMOGRAPHICS
tndoh_demographics <- readxl::read_excel(here::here("data_tndoh/Public-Dataset-RaceEthSex.xlsx"))
names(tndoh_demographics) <-  tolower(names(tndoh_demographics))
tndoh_demographics$cat_detail <- gsub("/ ", "/", tndoh_demographics$cat_detail)
tndoh_demographics <-
  tndoh_demographics %>%  
  filter(category == "RACE",
         date == max(date)) %>%  
  rename(race = cat_detail,
         cases_byrace = cat_totalcases) %>%  
  mutate(date = lubridate::ymd(date),
         race = factor(race,
                       levels = c("White", "Black or African American", "Asian",  "Other/Multiracial", "Pending",  "American Indian or Alaska Native", "Native Hawaiian or Other Pacific Islander"),
                       labels = c("White", "Black", "Asian",  "Other", "Pending",  "Am.Indian/Alaskan", "Hawaiian/Pacific")))

# US AND STATE DATA ####
# us_data <- tigris::states(cb = TRUE)
# us_data <-   sf::st_as_sf(us_data)
# names(us_data) <- tolower(names(us_data))
# us_data <-
#   us_data %>%  
#   filter(stusps %in% state.abb)
# saveRDS(us_data, here::here("data_rds/us_data.rds"))
# us_data <-  readRDS(here::here("data_rds/us_data.rds"))
# 
# 
# 
# # state_pop <- readxl::read_excel(here::here("state_populations.xlsx"))
# # state_pop$region <- gsub("\\.", "", state_pop$region)
# # state_pop <-
# #   state_pop %>% 
# #   pivot_longer(cols = -c(region, census_2010),
# #                names_to = "pop_year",
# #                names_prefix = "pop_",
# #                values_to = "pop_estimate")
# # saveRDS(state_pop, here::here("data_rds/state_pop.rds"))
# 
# state_pop <- readRDS( here::here("data_rds/state_pop.rds"))
# state_pop_2019 <-
#   state_pop %>%  
#   filter(pop_year == 2019) %>%  
#   filter(region %in% state.name)


# US DATA FROM THE COVIDTRACKING PROJECT ####
# COVIDTRACKING PROJECT QUIT TRACKING DATA ON MARCH 7, 2021; 
# DATA SOURCE CHANGED TO CDD; TEST POSITIVITY RATES NOT AVAILABLE THERE
# covus_currentdata <- readr::read_csv(here::here("data_covidtracking/daily.csv"))
# covus_currentdata <-
#   covus_currentdata %>%  
#   filter(state %in% state.abb,
#          !state %in% c("AK", "HI") ) %>%  
#   mutate(date = lubridate::ymd(date))

# covus_currentdata <- readRDS(here::here("data_covidtracking/covus_current.rds"))
# covus_currentdata <-
#   covus_currentdata %>%  
#   left_join(us_data,
#             by = c("state" = "stusps")) %>%  
#     left_join(state_pop_2019, by = c("name" = "region")) %>%  
#   mutate(test_density_pct = total/pop_estimate*100,
#          case_density_pct = positive/pop_estimate*100)  
# 
# covus_currentdata <-
#   covus_currentdata %>% 
#   group_by(state) %>%  
#   arrange(date) %>%  
#   mutate(new_cases_state_7davg = round(RcppRoll::roll_meanr(positiveIncrease, 7), 1),
#          new_deaths_state_7davg = round(RcppRoll::roll_meanr(deathIncrease, 7), 1),
#          new_deaths_state_7davg_100k = new_deaths_state_7davg/pop_estimate*100000,
#          new_testpos_pct   = round(positiveIncrease/totalTestResultsIncrease*100, 1),
#          new_testpos_pct_7davg   = round(RcppRoll::roll_meanr(positiveIncrease, 7)/RcppRoll::roll_meanr(totalTestResultsIncrease, 7)*100, 1),
#          new_cases_state_7davg_100k = new_cases_state_7davg/pop_estimate*100000
#   ) %>%  
#   ungroup()
# 
# covus_currentdata <-
#   covus_currentdata %>% 
#   group_by(date) %>%  
#   mutate(rank_newcases7davg = rank(-new_cases_state_7davg_100k)) %>%  
#   ungroup()
# 
# state_tbl <- 
#   tibble(state_name = state.name, 
#          state = state.abb)
# 
# covus_currentdata_2163 <-
#   covus_currentdata %>%  
#   mutate(geometry = st_transform(geometry, 2163)) %>%  
#     left_join(state_tbl, by = "state") %>%  
#   mutate(tt = sprintf("<tr><td><b><i>%s (%s)</i></b></td></tr>
#                      <tr><td>New Cases/100k (7d Avg): %s</td></tr>
#                      <tr><td>Test %% Positive (7d Avg): %1.1f%%</td></tr>
#                      <tr><td>State Rank (Cases/100k 7d Avg): %s</td></tr>
#                      <tr><td>Total Cases: %s</td></tr>
#                      <tr><td>Total Deaths: %s</td></tr>
#                      <tr><td>New Deaths/d/100k (7d Avg): %1.2f</td></tr>
#                      <tr><td>CFR: %1.2f%%</td></tr>
#                      <tr><td>Case Density: %1.1f%% of pop.</td></tr>
#                      <tr><td>Testing Density: %1.1f%% of pop.</td></tr>", 
#                      toupper(state_name),  paste0(month.name[lubridate::month(date)], " ", lubridate::mday(date)),
#                      comma(new_cases_state_7davg_100k, digits = 1), new_testpos_pct_7davg,  rank_newcases7davg, comma(positive, digits = 0),
#                      comma(death, digits = 0), new_deaths_state_7davg_100k, death/positive*100, case_density_pct, test_density_pct)) 

# US DATA FROM THE CDC ####
# NOTE: COVIDTRACKING PROJECT QUIT TRACKING DATA ON MARCH 7, 2021; 
# TEST POSITIVITY RATES NOT AVAILABLE THERE
# CDC DATA
cdc_currentdata <- readr::read_csv(here::here("data_cdc/cdc_state_info.csv"))
names(state.name) <- state.abb
cdc_currentdata <-
  cdc_currentdata %>% 
  filter(state %in% state.abb) %>% 
         # !state %in% c("AK", "HI") ) %>% 
  mutate(date = lubridate::mdy(submission_date),
         region = state.name[state])


cdc_currentdata <-
  cdc_currentdata %>%  
  # left_join(us_data,
  #           by = c("state" = "stusps")) %>%  
  left_join(state_pop_2019, by = "region") 


cdc_currentdata <-
  cdc_currentdata %>% 
  group_by(state) %>%  
  arrange(date) %>%  
  mutate(case_density_pct = tot_cases/pop_estimate*100,
         new_cases_state_7davg = round(RcppRoll::roll_meanr(new_case, 7), 1),
         new_deaths_state_7davg = round(RcppRoll::roll_meanr(new_death, 7), 1),
         new_deaths_state_7davg_100k = new_deaths_state_7davg/pop_estimate*100000,
         new_cases_state_7davg_100k = new_cases_state_7davg/pop_estimate*100000
  ) %>%  
  ungroup()

cdc_currentdata <-
  cdc_currentdata %>% 
  group_by(date) %>%  
  mutate(rank_newcases7davg = rank(-new_cases_state_7davg_100k)) %>%  
  ungroup()



state_tbl <-
  tibble(state_name = state.name,
         state = state.abb)


cdc_currentdata_2163 <-
  cdc_currentdata %>%  
  filter(date == max(date)) %>% 
  left_join(us_data, by = c("state" = "stusps"))  %>% 
  mutate(geometry = st_transform(geometry, 2163)) %>%  
  left_join(state_tbl, by = c("state", "state_name")) %>% 
  mutate(tt = sprintf("<tr><td><b><i>%s (%s)</i></b></td></tr><br>
                     <tr><td>New Cases/100k (7d Avg): %s</td></tr><br>
                     <tr><td>State Rank (Cases/100k 7d Avg): %s</td></tr><br>
                     <tr><td>Total Cases: %s</td></tr><br>
                     <tr><td>Total Deaths: %s</td></tr><br>
                     <tr><td>New Deaths/d/100k (7d Avg): %1.2f</td></tr><br>
                     <tr><td>CFR: %1.2f%%</td></tr><br>
                     <tr><td>Case Density: %1.1f%% of pop.</td></tr>",
                     toupper(state_name),  paste0(month.name[lubridate::month(date)], " ", lubridate::mday(date)),
                     comma(new_cases_state_7davg_100k, digits = 1), rank_newcases7davg, comma(tot_cases, digits = 0),
                     comma(tot_death, digits = 0), new_deaths_state_7davg_100k, tot_death/tot_cases*100, case_density_pct)) 
  # mutate(tt = sprintf("<tr><td><b><i>%s (%s)</i></b></td></tr>
  #                    <tr><td>New Cases/100k (7d Avg): %s</td></tr>
  #                    <tr><td>State Rank (Cases/100k 7d Avg): %s</td></tr>
  #                    <tr><td>Total Cases: %s</td></tr>
  #                    <tr><td>Total Deaths: %s</td></tr>
  #                    <tr><td>New Deaths/d/100k (7d Avg): %1.2f</td></tr>
  #                    <tr><td>CFR: %1.2f%%</td></tr>
  #                    <tr><td>Case Density: %1.1f%% of pop.</td></tr>",
  #                    toupper(state_name),  paste0(month.name[lubridate::month(date)], " ", lubridate::mday(date)),
  #                    comma(new_cases_state_7davg_100k, digits = 1), rank_newcases7davg, comma(tot_cases, digits = 0),
  #                    comma(tot_death, digits = 0), new_deaths_state_7davg_100k, tot_death/tot_cases*100, case_density_pct)) 

# county
tndoh_vax_county <-
  readxl::read_excel(here::here("data_tndoh/COVID_VACCINE_COUNTY_SUMMARY.XLSX"),
                       col_types = c("date", "text", rep("numeric", 6)))
    

names(tndoh_vax_county) <- tolower(names(tndoh_vax_county))
tndoh_vax_county <-
  tndoh_vax_county %>% 
  mutate(date = lubridate::ymd(date),
         county_lower = tolower(county)) %>% 
    left_join(TN_counties,
              by = "county_lower") %>% 
    rowwise() %>%  
    mutate(vaccine_count_pct = vaccine_count/population*100,
           recipient_count_pct = recipient_count/population*100) %>% 
    ungroup() %>% 
    filter(!is.na(statefp))

Testing Update- Interactive

interactive version- mouse/hover over bars for details.
Limited to Last 4 months for better detail

pl_new_tnsdc <-
  tndoh_daily %>% 
    filter(date >=  max(date) - 120) %>% 
  pivot_longer(cols = c(new_cases, new_hosp, new_deaths),
               names_to = "var", 
               values_to = "value") %>%  
  # filter(date > lubridate::ymd("2021-01-01")) %>% 
  ggplot(aes(x = date)) +
  geom_rect(aes(xmin = max(tndoh_daily$date)-7, xmax = max(tndoh_daily$date) + 1, ymin = -Inf, ymax = Inf),
            fill = "grey", alpha = 0.2) +
  ggiraph::geom_col_interactive(aes(y = value, fill = var, tooltip = tt_daily, data_id = date), alpha = 0.5,
                               position = "identity") +
  geom_line(aes(y = new_cases_7davg),color = "steelblue", size = 1) +
  labs(x=NULL,
       y="Results by day reported",
       subtitle = "Results by Day, New (Total)",
       caption = "Results within the past 7-10 days (shaded period)\nmay not be fully reported",
       fill = NULL) +
  scale_fill_manual(values = c("steelblue",
                               "red",
                               "orange"),
                    labels = c(paste0("Cases: ", comma(tndoh_daily$new_cases[tndoh_daily$date == max(tndoh_daily$date)], digits = 0),
                                      " (", comma(tndoh_daily$total_cases[tndoh_daily$date == max(tndoh_daily$date)], digits = 0), ")"),
                               paste0("Deaths: ", comma(tndoh_daily$new_deaths[tndoh_daily$date == max(tndoh_daily$date)], digits = 0),
                                      " (", comma(tndoh_daily$total_deaths[tndoh_daily$date == max(tndoh_daily$date)], digits = 0), ")"),
                               paste0("Hospitalizations: ", comma(tndoh_daily$new_hosp[tndoh_daily$date == max(tndoh_daily$date)], digits = 0),
                                      " (", comma(tndoh_daily$total_hosp[tndoh_daily$date == max(tndoh_daily$date)], digits = 0), ")")
                    )) +
  scale_y_continuous(expand = c(0,0),
                     sec.axis = sec_axis( trans=~./TN_pop2019*10^5, name="Cases/100k Population")) +
  scale_x_date(date_breaks = "1 month", date_labels = "%b" ) +
  theme(
    legend.position = c(.01, 1), 
    legend.justification = c(0, 1),
    legend.text = element_text(size = 10))

hover_css <- "background-color:orange;stroke:black;fill:orange;"
tooltip_css <- " background-color:rgba(255, 255, 255, 0.9);color:black;stroke:black;opacity:0.1;border:2px solid black;border-radius: 8px;font-family:Arial;font-size: 0.875em;"

ggiraph::girafe(ggobj = pl_new_tnsdc,
                width_svg = 8,
                options = list(
                  offx = 40, offy = 40,
                  opts_hover(css = hover_css),
                  opts_tooltip(css = tooltip_css)
                ) )

Today’s Overview

# PLOTS ####
pl_tests_new <-
  tndoh_daily %>% 
  # filter(date - 90) %>% 
  ggplot(aes(x=date)) +
  geom_col(aes(y = new_tests, fill = "Tests"), alpha = 0.5) +
  geom_col(aes(y = new_cases, fill = "Cases"), alpha = 0.7) +
  labs(x=NULL,
       fill = NULL,
       y="Tests by day reported",
       subtitle = "Tests by Day, New (Total)") +
  scale_fill_manual(values = c("Cases" = "steelblue",
                               "Tests" = "grey"),
                      labels = c(paste0("Cases:", comma(tndoh_daily$new_cases[tndoh_daily$date == max(tndoh_daily$date)], digits = 0),
                                        " (", comma(tndoh_daily$total_cases[tndoh_daily$date == max(tndoh_daily$date)], digits = 0), ")"),
                                 paste0("Tests:", comma(tndoh_daily$new_tests[tndoh_daily$date == max(tndoh_daily$date)], digits = 0), 
                                        " (", comma(tndoh_daily$total_tests[tndoh_daily$date == max(tndoh_daily$date)], digits = 0), ")"))) +
  scale_y_continuous(expand = c(0,0),
                     sec.axis = sec_axis( trans=~./TN_pop2019*10^5, name="per 100k Population")) +
  scale_x_date(date_breaks = "1 month", date_labels = "%b" ) +
  annotate(geom = "text", x=min(tndoh_daily$date)+1, y = 0.65*max(tndoh_daily$new_tests, na.rm = T), hjust = 0, 
           label = label_pos_pct, size = 3.5) +
  theme(#axis.text.x =  element_text(angle = 45, hjust = 1),
        legend.position = c(.01, 1), 
        legend.justification = c(0, 1),
        legend.text = element_text(size = 10))


# age ####

pl_age <-
  tndoh_age %>% 
  filter(date == max(date)) %>% 
  ggplot() +
# pl_age <-
#   pl_age +
  geom_blank(data = tndoh_age[tndoh_age$ar_totalcases == max(tndoh_age$ar_totalcases), ], 
             aes(x=age.f, y=ar_totalcases*1.1)) +
  geom_bar(aes(age.f, ar_totalcases), alpha = 0.5, fill = "steelblue", stat = "identity", ) +
  geom_bar(aes(age.f, ar_totaldeaths), alpha = 0.5, fill = "red", stat = "identity", ) +
  geom_text(aes(age.f, ar_totalcases*1.02, vjust =0, label = comma(ar_totalcases, digits = 0)),
            color = "steelblue", size = 3) +
    geom_text(aes(age.f, ar_totaldeaths*1.01, vjust =0,  label = comma(ar_totaldeaths, digits = 0)),
              color = "red", size = 3) +
    labs(x=NULL,
       y="Total Case Count",
       subtitle = "Total Cases by Age") +
  scale_y_continuous(expand = c(0,0)) +
  theme(axis.text.x =  element_text(angle = 45, hjust = 1)) 

# multi-axis plot
pl_age_death <-
  tndoh_age %>%  
    filter(date == max(date)) %>%  
  ggplot(aes(age.f, deaths_case_rate, group = date)) +
# pl_age_death <-
#   pl_age_death +
  geom_point(color = "red") +
  geom_line(color = "red")

# Dual axis 
ggp <- ggplot_build(pl_age)
ggp2 <- ggplot_build(pl_age_death)

ylim.prim <- ggp$layout$panel_params[[1]]$y.range   # in this example, precipitation
ylim.sec <- ggp2$layout$panel_params[[1]]$y.range    # in this example, temperatureb <- diff(ylim.prim)/diff(ylim.sec)
b <- diff(ylim.prim)/diff(ylim.sec)
a <- b*(ylim.prim[1] - ylim.sec[1])

pl_age_death_final <-
  tndoh_age %>% 
  filter(date == max(date)) %>%  
  ggplot(aes(x=age.f)) +
# pl_age_death_final <-
#   pl_age_death_final +
  geom_blank(data = tndoh_age[tndoh_age$ar_totalcases == max(tndoh_age$ar_totalcases), ], 
             aes(x=age.f, y=ar_totalcases*1.1)) +
  geom_bar(aes(age.f, ar_totalcases), alpha = 0.5, fill = "steelblue", stat = "identity", ) +
  geom_bar(aes(age.f, ar_totaldeaths), fill = "red", stat = "identity", ) +
  geom_line(aes(y = a + deaths_case_rate*b, group = date), color = "red", size = 1.5) +
  geom_point(aes(y = a + deaths_case_rate*b, group = date), color = "red", size = 2, shape = 21, fill = "white") +
  geom_text(aes(age.f, ar_totalcases*1.02, vjust =0, label = comma(ar_totalcases, digits = 0)),
            color = "steelblue", size = 3) +
    geom_text(aes(age.f, ar_totaldeaths*1.01, vjust =0,  label = comma(ar_totaldeaths, digits = 0)),
              color = "red", size = 3) +
    labs(x=NULL,
       y="Total Case Count",
       subtitle = "Total Cases by Age") +
  scale_y_continuous(expand = c(0,0),
                     sec.axis = sec_axis(~ (. - a)/b, name = "Crude Case Fatality (%)")) +
  theme(axis.text.x =  element_text(angle = 45, hjust = 1),
        axis.title.y.right = element_text(color = "red"),
        axis.text.y.right = element_text(color = "red")) 


race_total <- sum(tndoh_demographics$cases_byrace[tndoh_demographics$race != "Pending"])

# AGE HEATMAP
max_newcaserate_byage <- max(tndoh_age$ar_new_cases_100k_7davg, na.rm = T)
ar_newcase_heatmap_100k <-
  tndoh_age %>%  
  filter(date >= lubridate::ymd("2020-04-01")) %>%  
  ggplot(aes(x = date,
             y = age.f)) +
# ar_newcase_heatmap_100k <-
#   ar_newcase_heatmap_100k +
  geom_tile(aes(fill = ar_new_cases_100k_7davg,
                color = ar_new_cases_100k_7davg)) +
  viridis::scale_fill_viridis(labels = scales::comma,
                              limits = c(0, max_newcaserate_byage), oob = scales::squish,
                              option = "plasma", na.value="black") +
  viridis::scale_color_viridis(labels = scales::comma,
                               limits = c(0, max_newcaserate_byage), oob = scales::squish,
                               option = "plasma", na.value="black") +
  labs(fill = "New Cases/100k/d\n(7dAvg)",
       color  = "New Cases/100k/d\n(7dAvg)",
       x = NULL, y = NULL,
       subtitle = "New Case Rate, by Age (population adjusted)") +
  scale_x_date(date_breaks = "1 month", date_labels = "%b", expand = c(0, 0)) +
  theme_minimal() +
  theme(axis.text.y = element_text(size = 12),
        panel.grid = element_blank(),
        legend.position = "bottom",
        legend.background = element_blank(), 
        panel.border = element_blank(),
        legend.title = element_text(size = 10),
        legend.text =  element_text(size = 8))

# fig 1 ####  
pl_TOP <- cowplot::plot_grid(pl_new_tnsdc, pl_tests_new, 
                             pl_age_death_final, ar_newcase_heatmap_100k,
                             align = "h", axis = "bt",
                             nrow = 2)
cowplot::plot_grid(title_tnsdc, pl_TOP, nrow = 2, rel_heights = c(.1, 1))

# ggsave(here::here("figures/pl_tnmap.jpg"), dpi = 600, width = 12, height = 4, units = "in")

SUMMARY BY AGE In Tennessee, 1.83% of residents >80 and 0.74% of those age 71-80 have died of COVID-19, to date.

Age Group Population (n) Cases Cases
(% of Total) 		Deaths Deaths
(% of Total) 		Deaths/100k Deaths (% in Age) Given 1st Vaccine (% in Age)
0-10 821984 48740 6% 5 0% 0.61 0.00% NA%
11-20 855582 112276 13% 5 0% 0.58 0.00% 14.4%
21-30 942211 156758 18% 50 0% 5.31 0.01% 33.5%
31-40 887898 135251 16% 142 1% 15.99 0.02% 40.7%
41-50 850125 129221 15% 428 3% 50.35 0.05% 46.6%
51-60 893578 121777 14% 1125 9% 125.90 0.13% 55.3%
61-70 808164 86767 10% 2381 19% 294.62 0.29% 68.4%
71-80 516865 50023 6% 3818 30% 738.68 0.74% 75.8%
81+ 252767 25528 3% 4617 37% 1826.58 1.83% 67.3%
Population within each age group taken from US Census estimates for TN, 2019
Vaccination data last updated on 2021-07-01.

NOTES:

  • “new_test” numbers reported on “2020-07-13”, “2020-07-11”, and “2020-06-12” are not representative of the time trend (reclassification of how tests were reported, dumped on a a few dates) and are omitted
  • the reporting format changed to include “probable” cases on 6/12/2020. A total of 202 cases were classified as “probable” up to that date, and they were all included in the new case count on 6/12/2020.
  • Test reporting Note: Beginning on 6/12/2020 the number of all tests reported includes all tests, including repeat tests in the same person; an unclear number of tests were added on 6/12/2020 to account for this change in reporting.
  • COVID-19 data reported for June 29, 2020 reflect results from 2 days- due to an unplanned shutdown of the state surveillance system during the weekend.
  • Age strata populations are taken from the US Census estimate for TN, year 2019 estimate (the most recent available)
  • Test reporting note: on 4/19/2021 a large number of backlogged results were released and reflected in the 4/17-4/19 results; most cases were not active and reflect results from several weeks before this date.

Vaccinations

A bit of good news. Vaccination data is not updated every day- most recent data provided.

TN now provides detail on number of recipients and number of vaccinations.

COUNTY vaccination file from the state is empty today

latest_vax_date <- max(tndoh_vax_county$date, na.rm = T)
tndoh_vax_county %>% 
  filter(date == latest_vax_date,
         !is.na(statefp)) %>% 
    ggplot() +
    geom_sf(data = TN_counties) +
    geom_sf(aes(geometry = geometry,
                fill = recipient_count_pct),
            color = "grey30",
            inherit.aes = F) +
    geom_sf_text(aes(label=sprintf("%0.1f%%", recipient_count_pct), 
                     # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_deaths, pct_pos), 
                     geometry = geometry, 
                     color = cut(recipient_count_pct, 2)), 
                 size = 3) + 
    viridis::scale_fill_viridis(option = "C", direction = 1) +
    scale_color_manual(values = c("white", "black"), guide = F) +
    labs(fill = "%Pop.:",
         title = paste0("Vaccination Coverage (latest update on ", latest_vax_date, ")"),
         subtitle = sprintf("First Vaccinations: %s (%1.1f%% of TN Population)\nCompleted Vacc. : %s (%1.1f%% of Pop.)\nTotal Vaccinations: %s",
                            comma(max(tndoh_vax_state$recipient_count, na.rm = T), 0),
                            max(tndoh_vax_state$recipient_count_pct, na.rm = T),
                            comma(max(tndoh_vax_state$recip_fully_vacc, na.rm = T), 0),
                            max(tndoh_vax_state$recip_fully_vacc_pct, na.rm = T),
                            comma(max(tndoh_vax_state$vaccine_count, na.rm = T), 0))
         ) +
    theme_void() +
    theme(legend.position = "left")

Rolling 7-Day Averages & Top 10 Most Cases

note: The “New Hospitalizations” in this chart, as reported by the state, may not accurately reflect actual hospitalizations. Right y-axis is population adjusted (per 100k for TN). Ignore for “Test Positivity.”

Counties in the top 10 highest new case rate over the past 7 days; 7-day average trend-line for daily new cases

variable.labs <- c( "Cases", "Hospitalizations", "Deaths", "Tests", "Test Positivity (%)", "Vaccinations")
names(variable.labs) <- c( "new_cases", "new_hosp", "new_deaths", "new_tests", "pos_pct_today", "new_vaccine_count")

plot_state_trends <- 
  tndoh_daily %>%  
  mutate(pos_pct_today = pos_pct_today*100) %>%  
  # select(date, new_cases_7davg)
  # select(date, new_cases, new_deaths, new_tests, pos_pct_today) %>%  
  select(date, new_cases, new_hosp, new_deaths, new_tests, pos_pct_today, new_vaccine_count) %>% 
  pivot_longer(cols = -date,
               names_to = "variable",
               values_to = "value") %>%  
  group_by(variable) %>%  
  arrange(date) %>%  
  mutate(value_7davg = RcppRoll::roll_meanr(value, 7)) %>%  
  ggplot() +
# plot_state_trends <- 
#   plot_state_trends +
  geom_line(aes(date, value), size = 0.6, color = "grey") +
  geom_line(aes(date, value_7davg), size = 1, color = "red") +
  scale_y_continuous(sec.axis = sec_axis( trans=~./TN_pop2019*10^5, name="per 100k Population")) +
  labs(x=NULL,
       y="Rolling 7-day Average",
       title = "TN State Trends") +
  facet_wrap(~variable, scales = "free_y", nrow = 2,
             labeller = labeller(variable = variable.labs)) +
  theme(strip.text = element_text(size = 14))


# MIDTN counties plot ####
counties_top10 <-
  tndoh_countynew %>%   #select(new_cases_county_7davg)
  ungroup() %>%  
  filter(date == max(date, na.rm = T),
         name %in% TN_counties$name) %>%  # exclude the pending and out-of-state
  top_n(10, wt = new_cases_county_7davg) %>%  
  pull(name)

plot_top10_counties <- 
  tndoh_countynew %>%  
  filter(name %in% counties_top10) %>% 
  ggplot() +
# plot_top10_counties <- 
#   plot_top10_counties +
  # geom_bar(stat = "identity") +
  geom_col(aes(date, new_cases), fill = "steelblue", alpha = 0.5) +
  geom_col(aes(date, new_deaths), fill = "red") +
  geom_line(aes(date, new_cases_county_7davg), color = "steelblue", size = 1) +
  labs(x=NULL,
       y="New Results by day",
       fill = NULL,
       alpha = NULL,
       title = "Top 10 New Cases by Day (7d avg.)") +
  scale_y_continuous(expand = c(0,0)) +
  ylim(0,NA) +
  scale_x_date(date_breaks = "3 month", date_labels = "%b" ) +
  theme(#axis.text.x =  element_text(angle = 90, hjust = 1),
        legend.position = "none") +
  facet_wrap(~name, nrow = 2, 
             scales = "free_y")+
  theme(strip.text = element_text(size = 14))

cowplot::plot_grid(plot_state_trends,
                   plot_top10_counties,
                   nrow = 2, rel_heights = c(1.5,2))

Top 10: Largest Increase in Cases

Counties with the top 10 greatest weekly increase in New Cases, by % Change of 7-day average versus 7 days prior (limited to counties with average of >5 daily cases).

# MIDTN counties plot ####
counties_top10_increase <-
  tndoh_countynew %>%   
  ungroup() %>%  
  filter(date == max(date, na.rm = T),
       new_cases_county_7davg_1wkdelta > 0,
         new_cases_county_7davg >= 5,
         name %in% TN_counties$name) %>%  # exclude the pending and out-of-state
  top_n(n = 10, wt = new_cases_county_7davg_1wkdelta) %>%  
  pull(name)

plot_top10_increase <-
  tndoh_countynew %>%  
  filter(name %in% counties_top10_increase) %>% 
  ggplot() +
# plot_top10_increase + 
  # geom_bar(stat = "identity") +
  geom_col(aes(date, new_cases), fill = "steelblue", alpha = 0.5) +
  geom_col(aes(date, new_deaths), fill = "red") +
  geom_line(aes(date, new_cases_county_7davg), color = "steelblue", size = 1) +
  labs(x=NULL,
       y="New Results by day",
       fill = NULL,
       alpha = NULL,
       subtitle = "New Cases by Day: 10 Counties with the largest 1-week increase") +
  scale_y_continuous(expand = c(0,0)) +
  ylim(0,NA) +
  scale_x_date(date_breaks = "3 month", date_labels = "%b" ) +
  theme(#axis.text.x =  element_text(angle = 90, hjust = 1),
        legend.position = "none") +
  facet_wrap(~name, nrow = 2, 
             scales = "free_y")+
  theme(strip.text = element_text(size = 14))

plot_top10_increase

# ifelse(length(plot_top10_increase$data$county) == 0,
#        print("No increase in any counties"),
#        print(plot_top10_counties))

Best 10: Largest Decrease in Cases

Counties with the top 10 greatest weekly Decrease in New Cases, by % Change of 7-day average versus 7 days prior (limited to counties with average of >5 daily cases).

Only counties with a decline in case averages are shown (may be <10)

# MIDTN counties plot ####
counties_top10_decrease <-
  tndoh_countynew %>%
  ungroup() %>%
  filter(date == max(date, na.rm = T),
       new_cases_county_7davg_1wkdelta < 0,
           new_cases_county_7davg >= 5,
         name %in% TN_counties$name) %>%  # exclude the pending and out-of-state
  top_n(n = -10, wt = new_cases_county_7davg_1wkdelta) %>%
  pull(name)

# pl_best10 <-
#  tndoh_countynew %>%  
#   filter(name %in% counties_top10_decrease) %>% 
#   ggplot() +
# # pl_best10 <-
# #   pl_best10 +
#   # geom_bar(stat = "identity") +
#   geom_col(aes(date, new_cases), fill = "steelblue", alpha = 0.5) +
#   geom_col(aes(date, new_deaths), fill = "red") +
#   geom_line(aes(date, new_cases_county_7davg), color = "steelblue", size = 1) +
#   labs(x=NULL,
#        y="New Results by day",
#        fill = NULL,
#        alpha = NULL,
#        subtitle = "New Cases by Day: 10 Counties with the largest 1-week Decrease") +
#   scale_y_continuous(expand = c(0,0)) +
#   ylim(0,NA) +
#   scale_x_date(date_breaks = "3 month", date_labels = "%b" ) +
#   theme(#axis.text.x =  element_text(angle = 90, hjust = 1),
#         legend.position = "none",
#         strip.text = element_text(size = 14)) +
#   facet_wrap(~name, nrow = 2, 
#              scales = "free_y")
# pl_best10
# # cat("There were no counties with a decreased 7d average")


ifelse(length(counties_top10_decrease) >= 1,
       # IF >=1 with negative trend
       plot(pl_best10),
       # IF No counties with negative trend
       "There are no counties with a decreasing case trend")
## [1] "There are no counties with a decreasing case trend"

Interactive Tennessee Map

County color indicates the New Cases per day per county population (7 day average)

# make county data to list
county_plots <-
  tndoh_countynew %>% 
  filter(!is.na(date),
         date >= lubridate::ymd("2020-03-15")) %>%  
  group_by(name) %>%  
  nest() %>%  
  arrange(name) %>%  
  mutate(pl_new = map2(data, name,
                       ~ggplot(data = .x) +
                         geom_col(aes(date, new_cases, fill = "Cases"), alpha = 0.5) +
                         geom_col(aes(date, new_deaths, fill = "Deaths")) +
                         geom_line(aes(date, new_cases_county_7davg, fill = "Case 7-Day Avg"), color = "steelblue", size = 1) +
                         # geom_col(aes(date, new_hospitalized, fill = "Hospitalized")) +
                         labs(x=NULL,
                              y="New Results by day",
                              title = glue::glue("{.y} County"),
                              fill = NULL,
                              alpha = NULL,
                              subtitle = "New Results by Date Reported") + 
                         scale_fill_manual(values = c("Cases" = "steelblue",
                                                      "Case 7-Day Avg" = "steelblue",
                                                      "Deaths" = "red"),
                                           # "Hospitalized" = "orange"),
                                           labels = c(paste0("Total Cases: ", comma(.x$total_cases[.x$date == max(.x$date)], digits = 0)),
                                                      paste0("7d Avg: ", comma(.x$new_cases_county_7davg[.x$date == max(.x$date)], digits = 0)),
                                                      paste0("Total Deaths: ", comma(.x$total_deaths[.x$date == max(.x$date)], digits = 0))
                                                      # paste0("Hospitalized: ", .x$total_hospitalized[.x$date == max(.x$date)])
                                           )
                         ) +
                         scale_y_continuous(expand = c(0,0)) +
                         ylim(0,NA) +
                         scale_x_date(date_breaks = "1 month", date_labels = "%m/%d" ) +
                         theme(axis.text.x =  element_text(angle = 45, hjust = 1),
                               legend.position = "bottom")
  ))

county_data_byday_now <-
  left_join(TN_counties, 
            {tndoh_countynew %>%  
                filter(date == max(date, na.rm = T)) %>%   
                select(name, total_cases, new_cases_county_7davg, new_cases_per100kpop_7davg, new_cases_per100kpop_7davg_rank,
                       new_testpos_pct, new_testpos_pct_7davg, new_tests, new_pos_tests,
                       new_cases_rank, new_cases, new_cases_rank,
                       total_deaths, total_hospitalized)}) %>%  
  mutate(cases = replace_na(total_cases, 0),
         deaths = replace_na(total_deaths, 0)) %>%  
  left_join({county_plots %>%  select(name, pl_new)},
            by = "name")


lbls <-
  sprintf("<strong><i>%s Co.</i></strong><br/>
          Total Cases: %s<br/>
          New Cases/d/100k, 7dAvg (rank): %1.1f (%s/%s)<br/>
          New Cases Today (rank): %s (%s/%s)<br/>
          Test +%% Today: %1.1f%% (%s/%s)<br/>
          Test +%% (7dAvg): %1.1f%%<br/>
          Total Deaths: %s<br/>
          <i>Click to view time-course</i>", 
          toupper(county_data_byday_now$name), 
          comma(county_data_byday_now$total_cases, digits = 0), 
          county_data_byday_now$new_cases_per100kpop_7davg, county_data_byday_now$new_cases_per100kpop_7davg_rank, max(county_data_byday_now$new_cases_per100kpop_7davg_rank), 
          comma(county_data_byday_now$new_cases, digits = 0),county_data_byday_now$new_cases_rank, max(county_data_byday_now$new_cases_rank), 
          county_data_byday_now$new_testpos_pct, county_data_byday_now$new_pos_tests, county_data_byday_now$new_tests, 
          county_data_byday_now$new_testpos_pct_7davg,
          county_data_byday_now$total_deaths) %>%  
  lapply(htmltools::HTML)


# pal <- colorQuantile("viridis", domain = county_data_byday_now$new_cases_per100kpop_7davg, n = 10, reverse = FALSE)
pal <- colorNumeric("viridis", domain = county_data_byday_now$new_cases_per100kpop_7davg, reverse = FALSE)


leaflet() %>%  
  addProviderTiles(providers$Stamen.TonerLite) %>% 
  addPolygons(data = TN_counties, # borders of all counties
              color = "grey",
              fill = NA,
              weight = 1) %>% 
  addPolygons(data = county_data_byday_now,
              group = "name",
              fillColor = ~pal(county_data_byday_now$new_cases_per100kpop_7davg),
              fillOpacity = 0.8,
              weight = 1,
              highlight = highlightOptions(
                weight = 5,
                color = "red",
                bringToFront = TRUE),
              label = lbls,
              labelOptions = labelOptions(
                style = list("font-weight" = "normal", padding = "3px 8px"),
                textsize = "15px",
                direction = "auto")) %>% 
  addLegend(pal = pal, 
            values = county_data_byday_now$new_cases_per100kpop_7davg,
            title = "New Cases/100k (7d Avg)",
            opacity = 0.9,
            position = "bottomright") %>%  
  leafpop::addPopupGraphs(county_data_byday_now$pl_new,
                          group = "name",
                          width = 400, height = 200)

TN County Heatmap

Counties are ordered by highest new cases per day (current) based on the most recent 7 day average (per 100k population). Max limit of the color scale set at 200 to increase contrast- some local outbreaks exceed this limit (eg prison outbreaks in Bledsoe, Lake, Wayne Co.)

county_order <-
  tndoh_countynew %>%  
  filter(date == max(date)) %>%  
  mutate(county_bynewcases = reorder(county, new_cases_per100kpop_7davg)) %>%  
  select(county, county_bynewcases)

# no longer interactive due to file size issues.
TN_heatmap_100k <-
  tndoh_countynew %>%  
    filter(!county %in% c("Pending", "Out of State")) %>% 
  left_join(county_order, by = "county") %>%  
  mutate(new_cases = replace_na(new_cases, 0)) %>%  
  # mutate(lbls =
  #          sprintf("<tr><td><strong><i>%s Co. (%s)</i></strong></tr></td>
  #         <tr><td>Total Cases: %s</tr></td>
  #         <tr><td>New Cases/d/100k, 7dAvg (rank): %1.1f (%s/%s)</tr></td>
  #         <tr><td>New Cases Today (rank): %s (%s/%s)</tr></td>
  #        <tr><td> Test +%% Today: %1.1f%% (%s/%s)</tr></td>
  #         <tr><td>Test +%% (7dAvg): %1.1f%%</tr></td>
  #         <tr><td>Total Deaths: %s</tr></td>", 
  #         toupper(name), paste0(month.name[lubridate::month(date)], " ", lubridate::mday(date)),
  #         comma(total_cases, digits = 0), 
  #         new_cases_per100kpop_7davg, new_cases_per100kpop_7davg_rank, max(new_cases_per100kpop_7davg_rank, na.rm = T), 
  #         comma(new_cases, digits = 0),new_cases_rank, max(new_cases_rank, na.rm = T), 
  #         new_testpos_pct, new_pos_tests, new_tests, 
  #         new_testpos_pct_7davg,
  #         total_deaths) ) %>%  
  ggplot(aes(x = date,
             y = county_bynewcases,
             fill = new_cases_per100kpop_7davg)) +
# TN_heatmap_100k <-
#   TN_heatmap_100k +
  geom_tile(aes(fill = new_cases_per100kpop_7davg,
                color = new_cases_per100kpop_7davg)) +  
  # ggiraph::geom_tile_interactive(aes(tooltip = lbls,
  #                                    data_id = county,
  #                                    fill = new_cases_per100kpop_7davg,
  #                                    color = new_cases_per100kpop_7davg)) +
  viridis::scale_fill_viridis(labels = scales::comma,
                                 limits = c(0, n_ul), oob = scales::squish,
                              option = "plasma") +
  viridis::scale_color_viridis(labels = scales::comma,
                               limits = c(0, n_ul), oob = scales::squish,
                               option = "plasma") +
  labs(fill = "New Cases/100k\n(7dAvg.)",
       color = "New Cases/100k\n(7dAvg.)",
       x = NULL, y = NULL) +
  scale_x_date(date_breaks = "1 month", date_labels = "%b", expand = c(0, 0)) +
  theme_minimal() +
  theme(axis.text.y = element_text(size = 8),
        panel.grid = element_blank(),
        legend.position = "bottom",
        legend.background = element_blank(), 
        panel.border = element_blank(),
        legend.title = element_text(size = 14),
        legend.text =  element_text(size = 10))

# x_TN_heatmap_100k <- girafe(ggobj = TN_heatmap_100k,
#             options = list(
#               opts_hover(css = "stroke-width:0.5")))
TN_heatmap_100k

Notes:

  • 12/20/2020: Increased upper limit to 200 because 125 just wasn’t enough.

New cases and Hospitalizations by county

Note: assignment of new cases to county usually lags the total count. The State Health Dept. assigns cases to county of residence, rather than testing site. Individual county health departmens may report numbers differently, such as the number of cases diagnosed within the county, regardless of residence

# fig 2 State ####

pl_tnmap_tnsdc <-
  tndoh_countynew %>%  
  filter(!is.na(date)) %>%  
  filter(date == max(date)) %>% 
  mutate(total_cases = case_when(total_cases == 0 ~ NA_real_,
                               TRUE ~ total_cases)) %>%  
  ggplot() +
# pl_tnmap_tnsdc <-
#   pl_tnmap_tnsdc +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = as.numeric(total_cases)),
          inherit.aes = F) +
  geom_sf_text(aes(label= comma(total_cases, digits = 0), 
                   geometry = geometry, 
                   color = cut(total_cases, 2)), 
               size = 2.5) + 
  viridis::scale_fill_viridis(labels = scales::comma,
                                 # limits = c(0, 100), oob = scales::squish,
                              option = "plasma") +
    # viridis::scale_fill_viridis(option = "D", direction = -1) +
    scale_color_manual(values = c( "white", "black"), guide = F) +
  labs(fill = "Total Cases:")+
         # caption = "Visualization by @DrJMLuther\nData Source:@TNDeptofHealth via TN State Data Center") +
  theme_void() +
  theme(legend.position = "left")  
# ggsave(here::here("figures/pl_tnmap.jpg"), dpi = 600, width = 12, height = 4, units = "in")

pl_tnmap_7dDelta <-
  tndoh_countynew %>%  #select(new_cases_county_7davg_1wkdelta)
  filter(!is.na(date)) %>%  
  filter(date == max(date)) %>% 
  ggplot() +
# pl_tnmap_7dDelta <-
#   pl_tnmap_7dDelta +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = new_cases_county_7davg_1wkdelta*100),
          inherit.aes = F) +
  scale_fill_gradient2(low = "green", mid = "white", high = "red") +
  geom_sf_text(aes(label=sprintf("%1.1f",new_cases_county_7davg_1wkdelta*100), geometry = geometry),
               # color = cut(new_cases_county_7davg_1wkdelta, 2)),
               size = 3) +
  # viridis::scale_fill_viridis(option = "D", direction = -1) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "1-Week\nChange %:")+
  # caption = "Visualization by @DrJMLuther\nData Source:@TNDeptofHealth via TN State Data Center") +
  theme_void() +
  theme(legend.position = "left")  

pl_tnmap_newhosp <-
  tndoh_countynew %>%  
  filter(!is.na(date)) %>%  
  filter(date == max(date)) %>%  
  ggplot() +
# pl_tnmap_newhosp <-
#   pl_tnmap_newhosp +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = as.numeric(new_hospitalized)),
          inherit.aes = F) +
  geom_sf_text(aes(label=new_hospitalized, geometry = geometry), 
               color = "black",
               # color = cut(new_hospitalized, 2)), 
               size = 3) + 
  scale_fill_gradient(low = "white", high = "red", 
                      limits = c(0, max(tndoh_countynew$new_hospitalized [tndoh_countynew$date == max(tndoh_countynew$date, na.rm = T)], na.rm = T)))  +
  scale_color_manual(values = c( "white", "black"), guide = F) +
  labs(fill = "New\nHospitalizations:") +
  theme_void() +
  theme(legend.position = "left")  

pl_tnmap_new_tnsdc <-
  tndoh_countynew %>%  
  filter(date == max(date)) %>%  
  ggplot() +
# pl_tnmap_new_tnsdc <-
#   pl_tnmap_new_tnsdc +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = as.numeric(new_cases)),
          inherit.aes = F) +
  geom_sf_text(aes(label=new_cases, geometry = geometry), 
               color = "black",
               # color = cut(new_cases, 2)), 
               size = 3) + 
  scale_fill_gradient(low = "white", high = "red", 
                      limits = c(0, max(tndoh_countynew$new_cases[tndoh_countynew$date == max(tndoh_countynew$date, na.rm = T)], na.rm = T)))  +
  labs(fill = "New Cases:") +
       # caption = "Visualization by @DrJMLuther\nData Source:@TNDeptofHealth via New York Times") +
  theme_void() +
  theme(legend.position = "left")  
# ggsave(here::here("figures/pl_tnmap_new.jpg"), dpi = 600, width = 12, height = 4, units = "in")


pl_FINAL_tnsdc <-
  cowplot::plot_grid(title_tnsdc,
                   pl_tnmap_new_tnsdc,
                   pl_tnmap_7dDelta,
                   pl_tnmap_newhosp,
                   # pl_tnmap_newtests_tnsdc,
                   pl_tnmap_tnsdc,
                   nrow = 5,
                   rel_heights = c(.1,1,1,1, 1),
                   align = "y")
  unclassified <- 
    tndoh_countynew %>%  
    filter(date == max(date),
         county_lower %in% c("pending", "out of state"))

  
    # pl_FINAL_tnsdc
  pl_FINAL_tnsdc +
    annotation_custom(gridExtra::tableGrob(
      tibble(
        Category = unclassified$county,
        Cases = unclassified$new_cases), 
      rows=NULL, cols = NULL,
      theme = gridExtra::ttheme_minimal(base_size = 10)),
      xmin=.85, xmax=.95, ymin=0.66, ymax=.9) +
    annotation_custom(gridExtra::tableGrob(
      tibble(
        Category = unclassified$county,
        Cases = unclassified$new_hospitalized), 
      rows=NULL, cols = NULL,
      theme = gridExtra::ttheme_minimal(base_size = 10)),
      xmin=.85, xmax=.95, ymin=0.31, ymax=.275) +
    annotation_custom(gridExtra::tableGrob(
      tibble(
        Category = unclassified$county,
        Cases = unclassified$total_cases), 
      rows=NULL, cols = NULL,
      theme = gridExtra::ttheme_minimal(base_size = 10)),
      xmin=.85, xmax=.95, ymin=0.01, ymax=.125)

Deaths and Crude Case Fatality ratio

The estimate of case fatality is presented as the crude deaths/cases ratio, and is not adjusted for age. This could markedly affect the rate in different counties. This is also affected by the ability to test and identify as many cases as possible, which may differ by county.

pl_deaths <-
  tndoh_countynew %>%    
  filter(date == max(date, na.rm = T)) %>%  
  ggplot() +
# pl_deaths <-
#   pl_deaths +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = total_deaths),
          color = "grey30",
          inherit.aes = F) +
  geom_sf_text(aes(label=sprintf("%0.0f", comma(total_deaths, digits = 0)), 
  # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_deaths, pct_pos), 
                   geometry = geometry, 
                   color = cut(total_deaths, 2)), 
               size = 3) + 
  viridis::scale_fill_viridis(option = "C", direction = 1) +
  scale_color_manual(values = c("white", "black"), guide = F) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "Deaths:") +
  theme_void() +
  theme(legend.position = "left") 


pl_cfr <-
  tndoh_countynew %>%    
  filter(date == max(date, na.rm = T)) %>%  
  mutate(cfr = total_deaths/total_cases) %>%  
  ggplot() +
# pl_cfr <-
#   pl_cfr +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = cfr*100),
          color = "grey30",
          inherit.aes = F) +
  geom_sf_text(aes(label=sprintf("%0.1f %%", cfr*100), 
  # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_cases, pct_pos), 
                   geometry = geometry, 
                   color = cut(cfr, 2)), 
               size = 3) + 
  viridis::scale_fill_viridis(option = "C", direction = 1) +
  scale_color_manual(values = c("white", "black"), guide = F) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "Case Fatality %:") +
  theme_void() +
  theme(legend.position = "left") 




cowplot::plot_grid(title_tnsdc, pl_deaths, pl_tnmap_tnsdc, pl_cfr,
                   nrow = 4, rel_heights = c(.1,1,1, 1))

# ggsave(here::here("figures/pl_deaths_map.jpg"), dpi = 600, width = 12, height = 6, units = "in")

County Population adjusted rates

The estimated Tennessee county populations are taken from the US Census (2019 estimate). The estimate of case fatality is presented as the crude deaths/cases ratio, and is not adjusted for age. This could markedly affect the rate in different counties

pl_deaths <-
  tndoh_countynew %>%    
  filter(date == max(date, na.rm = T)) %>%  
  ggplot() +
# pl_deaths <-
#   pl_deaths +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = total_deaths_per100kpop),
          color = "grey30",
          inherit.aes = F) +
  geom_sf_text(aes(label=sprintf("%0.0f", total_deaths_per100kpop), 
                   # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_deaths, pct_pos), 
                   geometry = geometry, 
                   color = cut(total_deaths_per100kpop, 2)), 
               size = 3) + 
  viridis::scale_fill_viridis(option = "C", direction = 1) +
  scale_color_manual(values = c("white", "black"), guide = F) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "Deaths per\n100K Pop.:") +
  theme_void() +
  theme(legend.position = "left") 

pl_cases <-
  tndoh_countynew %>%    
  filter(date == max(date, na.rm = T)) %>%  
  ggplot() +
# pl_cases <-
#   pl_cases +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = total_cases_per100kpop),
          color = "grey30",
          inherit.aes = F) +
  geom_sf_text(aes(label=sprintf("%0.1f%%", total_cases_per100kpop/1000), 
                   # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_cases, pct_pos), 
                   geometry = geometry, 
                   color = cut(total_cases_per100kpop, 2)), 
               size = 3) + 
  viridis::scale_fill_viridis(option = "C", direction = 1) +
  scale_color_manual(values = c("white", "black"), guide = F) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "Population % Affected\n(Cases %Pop.):") +
  theme_void() +
  theme(legend.position = "left") 

pl_newcases <-
  tndoh_countynew %>%  #select(new_cases_per100kpop_7davg)   
  filter(date == max(date, na.rm = T)) %>%  
  ggplot() +
# pl_newcases <-
#   pl_newcases +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = new_cases_per100kpop_7davg),
          color = "grey30",
          inherit.aes = F) +
  geom_sf_text(aes(label=sprintf("%0.0f", new_cases_per100kpop_7davg), 
                   # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_cases, pct_pos), 
                   geometry = geometry, 
                   color = cut(new_cases_per100kpop_7davg, 2)), 
               size = 3) + 
  viridis::scale_fill_viridis(option = "C", direction = 1) +
  scale_color_manual(values = c("white", "black"), guide = F) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "Daily Cases per\n100K Pop.\n(7d Avg):") +
  theme_void() +
  theme(legend.position = "left") 

pl_cfr <-
  tndoh_countynew %>%  
    mutate(total_tests_pctpop = total_tests/population*100000) %>%  
    # select(total_tests_per100kpop)   
  filter(date == max(date, na.rm = T)) %>%  
  ggplot() +
# pl_cfr <-
#   pl_cfr +
  geom_sf(data = TN_counties) +
  geom_sf(aes(geometry = geometry,
              fill = total_tests_pctpop),
          color = "grey30",
          inherit.aes = F) +
  geom_sf_text(aes(label=sprintf("%0.1f ", total_tests_pctpop/100), 
                   # geom_sf_text(aes(label=sprintf("%0.0f\n%0.0f%%", total_cases, pct_pos), 
                   geometry = geometry, 
                   color = cut(total_tests_pctpop, 2)), 
               size = 3) + 
  viridis::scale_fill_viridis(option = "C", direction = 1) +
  scale_color_manual(values = c("white", "black"), guide = F) +
  # scale_color_manual(values = c( "black", "white"), guide = F) +
  labs(fill = "Tests per\n100k Population:") +
  theme_void() +
  theme(legend.position = "left") 

cowplot::plot_grid(  title_tnsdc, pl_newcases,pl_deaths, pl_cases,pl_cfr,
                   nrow = 5, rel_heights = c(.1,1,1,1, 1))

# ggsave(here::here("figures/pl_deaths_map.jpg"), dpi = 600, width = 12, height = 6, units = "in")

US Data from the CDC

Case data for these maps are taken from the CDC. and the 2019 state population estimates are taken from the United States Census Bureau.

Formerly I used the COVID Tracking project, but this site is no longer updating data as of March 7, 2021.

New Cases

mouse over state for details

max_newcaserate <- max(cdc_currentdata_2163$new_cases_state_7davg_100k[cdc_currentdata_2163$date == max(cdc_currentdata_2163$date)], na.rm = T)
                                                  
pl_us_cases_100k <-
  cdc_currentdata_2163 %>%  
  filter(date == max(date)) %>%  
  ggplot(aes(geometry = geometry,
             fill = new_cases_state_7davg_100k)) +
# pl_us_cases_100k <-
#   pl_us_cases_100k +
  ggiraph::geom_sf_interactive(aes(tooltip = tt,
                                   data_id = state), size = .2) +
  viridis::scale_fill_viridis(labels = scales::comma,
                               limits = c(0, max_newcaserate), oob = scales::squish,
                              option = "plasma") +
  labs(fill = "Cases/100k/d\n(7dAvg)") +
  theme_void() 

ggiraph::girafe(code = print( pl_us_cases_100k),
                options = list(
                  opts_hover(css = "stroke:orange;stroke-width:2;")))

US Heatmap

max_newcaserate_alldates <- max(cdc_currentdata$new_cases_state_7davg_100k, na.rm = T)

state_order <- 
  cdc_currentdata_2163 %>%  
  filter(date == max(date)) %>%  
  mutate(state_bynewcases = reorder(state, new_cases_state_7davg_100k)) %>%  
  select(state, state_bynewcases)

us_heatmap_100k <-
  cdc_currentdata %>%  
  left_join(state_order, by = "state") %>%  
  filter(date >= lubridate::ymd("2020-04-01")) %>%  
  ggplot(aes(x = date,
             y = state_bynewcases)) +
# us_heatmap_100k <-
#   us_heatmap_100k +
  geom_tile(aes(fill = new_cases_state_7davg_100k,
                color = new_cases_state_7davg_100k)) +
  # ggiraph::geom_tile_interactive(aes(tooltip = tt,
  #                                    data_id = state,
  #                                    fill = new_cases_state_7davg_100k,
  #                                    color = new_cases_state_7davg_100k)) +
  viridis::scale_fill_viridis(labels = scales::comma,
                               limits = c(0, max_newcaserate_alldates), oob = scales::squish,
                              option = "plasma") +
    viridis::scale_color_viridis(labels = scales::comma,
                               limits = c(0, max_newcaserate_alldates), oob = scales::squish,
                              option = "plasma") +
  labs(fill = "Cases/100k/d\n(7dAvg)",
       color  = "Cases/100k/d\n(7dAvg)",
       x = NULL, y = NULL) +
  scale_x_date(date_breaks = "1 month", date_labels = "%b", expand = c(0, 0)) +
  theme_minimal() +
  theme(axis.text.y = element_text(size = 12),
        panel.grid = element_blank(),
        legend.position = "bottom",
        legend.background = element_blank(), 
        panel.border = element_blank(),
        legend.title = element_text(size = 14),
        legend.text =  element_text(size = 12))
us_heatmap_100k

You made it to the end, so here’s a Dad joke:

I used to hate facial hair…but then it grew on me.

Session Info

Hmisc::markupSpecs$html$session()
 R version 4.1.0 (2021-05-18)
 Platform: x86_64-w64-mingw32/x64 (64-bit)
 Running under: Windows 10 x64 (build 17763)
 
 Matrix products: default
 
 attached base packages:
 [1] stats     graphics  grDevices utils     datasets  methods   base     
 
 other attached packages:
  [1] urbnmapr_0.0.0.9002 ggiraph_0.7.10      formattable_0.2.1  
  [4] leafpop_0.1.0       leaflet_2.0.4.1     sf_1.0-1           
  [7] forcats_0.5.1       stringr_1.4.0       dplyr_1.0.7        
 [10] purrr_0.3.4         readr_1.4.0         tidyr_1.1.3        
 [13] tibble_3.1.2        ggplot2_3.3.5       tidyverse_1.3.1
To cite R in publication use:

R Core Team (2021). R: A Language and Environment for Statistical Computing. R Foundation for Statistical Computing, Vienna, Austria. https://www.R-project.org/.